Essa é uma revisão anterior do documento!
Fisica Computacional
Professor: Marcio Argollo
Contato: Sala A1-02, tel. 2629-5773
Referências básicas:
- Manuais de Fortran 77:
- Manuais de Fortran 90
- Fortran 90 para programadores de Fortran 77
- Manuais de C/C++
- Livros-texto
- An Introduction to Computer Simulation Methods: Applications to Physical Systems, H. Gould e J. Tobochnik.
Para uma revisão de comandos básicos em Unix/Linux, utilização do gnuplot e programação em C veja o tutorial apresentado na Semana da Física 2009 "Ferramentas de Computação em Física".
Listas de Exercícios
As listas deverão ser entregues em papel e por email, compactadas no formato tar.gz e com o nome listaN_seunome.tar.gz
Lista 1) Exercicios de revisão de programação - 19/08/2010
- Calcular as somas
e 
com i sendo uma variável ponto flutuante de precisão simples. Plote as diferencas 
vs N e 
vs N para 
onde 
- Escrever um programa que recebe como entrada um número inteiro em representação decimal e retorna sua expansão binária. Para tal, crie uma função que verifica se o i-ésimo bit de uma palavra (variável do programa) vale zero ou um.
- Considere um número de ponto flutuante de N bits com E bits reservados para o expoente e M=N-E-1 para a mantissa. Represente TODOS os números obtidos nesta representação em um grafico de barras no gnuplot. Para cada número x, armazene um par (x,1) e plote o arquivo de dados com
with impulses. Reveja aqui a definição de números de ponto flutuante.
- Imagine o empilhamento de dominós de massa unitária e comprimento 2. O primeiro dominó se encontra com o centro de massa na origem, e a cada passo n deslocamos os dominós ja empilhados de
para a direita e adicionamos outro dominó na parte inferior da pilha, também com seu centro de massa na origem. Escreva um programa que diz, para um dado 
, quantos dominós podem ser empilhados até que o centro de massa ultrapasse a largura do primeiro dominó e a pilha tombe. (Veja uma discussão do problema aqui).
Lista 2) Derivada numérica e integração de equações diferenciais ordinárias
O método de Euler - 31/08/2010
Sendo a derivada de uma função 
definida por 
e
, o método de Euler (com difereças posteriores) corresponde ao truncamento desta expansão em ordem h. Obtemos assim que 
.
- Mostre que o método de Euler com diferenças centradas é exato em
, onde h é o passo de discretização da derivada. Nesse método escrevemos 
- Calcule a derivada numérica de
no ponto x=3 usando o método de Euler simples com passo de discretização h entre 
e e faça um gráfico mostrando a diferença relativa entre a derivada numérica e o valor exato da derivada para ponto flutuante com precisão simples e dupla. Entre quais valores de h a aproximação é aceitável? - Mostre que o método de Euler é instável para qualquer passo h na solução do problema de crescimento exponencial
- Integre a equação diferencial
com 
com o método de Euler simples para diferentes valores de passo h=0.1, h = 0.05, h = 0.01, h = 0.005 e h = 0.001 e compare cada solução em t=1,2,3,4 e 5 com a solução exata 
. Faça um gráfico com a solução para h=0.05 entre t=0 e t=5 e a solução exata. - Resolva o problema acima com o método de Euler com diferenças centradas e os mesmos valores de h. Para obter o primeiro passo use o método de Euler simples com passo de discretização dez vezes menor que o usado no resto do intervalo.
- Referências
Integração numérica de equações diferenciais ordinárias de sistemas conservativos: métodos simpléticos - 01/09/2010
A aproximação para a derivada, como sugerida no método de Euler, introduz erros sistemáticos que acarretam na destruição de certos invariantes da dinâmica original, como a energia total no movimento do pêndulo simples. Para reduzir os efeitos da aproximação numérica à derivada, pode-se recorrer a métodos com aproximações de ordem superior (o método de Euler é de primeira ordem no passo de discretização h).
Métodos de Runge-Kutta
Ver: http://mathworld.wolfram.com/Runge-KuttaMethod.html
- Resolva numericamente o problema do pêndulo físico utilizando o método de Runge-Kutta de quarta ordem com passo
e ilustre o espaço de fase do sistema. Considere 
metros, onde é o comprimento da haste de massa desprezível à qual se prende uma massa m. Use a condição inicial 
(posição vertical para baixo) e velocidades iniciais 
e 
, onde 
é a velocidade crítica que separa os movimentos de oscilação e rotação (Ver o artigo "Comportamento crítico no pêndulo simples" do Paulo Murilo para um estudo desta transição). Considere a evolução por dois períodos do movimento. Ilustre, no mesmo gráfico, a solução do problema na aproximação de pequenos ângulos com e 
. - Faça um gráfico contendo os primeiros 30 segundos da evolução temporal dos valores analítico e numérico da energia total do pêndulo. No caso numérico, utilize o método de Euler e o método de Runge-Kutta de segunda e quarta ordem com passos
e 
.
Métodos simpléticos
Ver Reversible multiple time scale molecular dynamics (pdf aqui)
Estes métodos procuram criar esquemas de discretização da evolução Hamiltoniana de um sistema conservativo, de modo que o volume de pontos representativos no espaço de fase seja preservado durante a evolução temporal. Estes esquemas envolvem fatorização do operador Liouvilleano.
- Use o método de Euler-Cromer e o método de Euler com passo e para resolver numericamente o problema do pêndulo simples na aproximação de pequenas oscilações. Ulitize a condição inicial
e 
. Faça um gráfico contendo a evolução temporal da posição angular por aproximadamente 10 períodos do movimento para a solução analítica e numérica com os métodos de Euler e Euler-Cromer. Faça outro gráfico ilustrando a evolução temporal da energia no mesmo alcance temporal para as 3 soluções.
Lista 3: Caos em sistemas dissipativos: o pêndulo forçado-amortecido -14/09/2010
Considere o pêndulo simples em um meio viscoso, que oferece resistência ao seu movimento com uma força proporcional à sua velocidade 
. A equação de movimento pode ser escrita como
.
- Faça o gráfico da posição como função do tempo para
e 
com posição angular inicial 
e 
. Note que para valores grandes de q devemos recuperar o movimento sem dissipação (conservativo) e, portanto, devemos utilizar um método simplético, como Verlet ou Euler-Cromer, para integrar a equação de movimento.
Considere agora uma força periódica impulsionando o pêndulo com força 
. A equação de movimento pode ser escrita como
.
- Construa o gráfico
com 
e 
e 
, respectivamente. Observe a mudança de comportamento ao variarmos a amplitude A da força externa. Simule o pêndulo por 100 ciclos da força externa e despreze os primeiros 10 ciclos (quanto vale o período 
da força externa neste exemplo?). - Construa o espaço de fases para cada uma das situações descritas acima com a evolução do ponto representativo
- Faça um gráfico “estroboscópico”, chamado seção de Poincaré, imprimindo os pontos no espaço de fase correspondentes a intervalos de tempo múltiplos do período da força externa
. - Obtenha uma série de pontos
onde 
corresponde ao valor da velocidade angular quando 
vale zero pela n-ésima vez durante a evolução temporal do pêndulo. Despreze os valores obtidos para tempos menores que 10 ciclos da força externa. Construa agora outra seção de Poincaré em um gráfico 
para os mesmos parâmetros do problema anterior.
e 
, obtemos comportamento caótico para alguns valores da amplitude da força externa A. Verifique!
| Amplitude da força | Comportamento |
|---|---|
| A < 1.085 | periódico |
| 1.085 < A < 1.11 | caótico |
| 1.11 < A < 1.14 | periódico |
| 1.14 < A < 1.22 | caótico |
| A ~ 1.22 | periódico |
| 1.22 < A < 1.28 | caótico |
| 1.28 < A < 1.475 | periódico |
| 1.475 < A < 1.485 | caótico |
| 1.485 < A < 1.493 | periódico |
| 1.493 < A < 1.495 | caótico |
| 1.495 < A < 1.497 | periódico |
Lista 4: Caos no mapa logístico - 23/09/2010
Leia aqui sobre o problema.
Considere o mapa logístico 
- Faça os gráficos (escala lin-log) da evolução temporal de x a partir de um valor inicial arbitrário para 
e 
.
- Fazer “gráficos de escada” para os mesmos valores de parâmetro e condição incial x=0.6. Como critério de parada use a condição 
.
- Fazer “gráficos de escada” para o caso 
e condição inicial x=0.3 e para a segunda iterada f(f(x)) com mesma condição inicial.
- Construir o diagrama de bifurcação do mapa logístico na região 
.
- Obter os expoentes de Lyapunov 
do mapa logístico na região .
O mapa logístico com 
é ergódico, o que significa que podemos substituir médias temporais por médias sobre configurações de um ensemble e, nesse caso,
onde 
é a fração de vezes que a órbita se encontra entre x e x+dx.
- Calcule, subdividindo o intervalo [0,1] em 100 subintervalos, o histograma para as primeiras 100000 iterações do mapa logístico com e encontre . Esse caso particular admite a solução analítica 
. Compare os resultados mostrando-os em um mesmo gráfico.
- Obter o diagrama de bifurcação para os mapas
(a) 
(b) 
ambos na região 
Lista 5: Análise de séries temporais - 19/10/2010
Voce encontra neste arquivo a sequência dos primeiros 10 mil dígitos de 
.
- Calcular a frequência de cada dígito i na sequência.
- Calcular a frequência de todos os pares de dígitos (i,j) na sequência.
